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Abstract 

We develop a multi-factor stochastic volatility Libor model with dis- 
placement, where each individual forward Libor is driven by its own 
square-root stochastic volatility process. The main advantage of this ap- 
proach is that, maturity- wise, each square-root process can be calibrated 
to the corresponding cap(let)vola-strike panel at the market. However, 
since even after freezing the Libors in the drift of this model, the Libor 
dynamics are not affine, new affine approximations have to be developed in 
order to obtain Fourier based (approximate) pricing procedures for caps 
and swaptions. As a result, we end up with a Libor modeling package 
that allows for efficient calibration to a complete system of cap/swaption 
market quotes that performs well even in crises times, where structural 
breaks in vola-strike-maturity panels are typically observed. 
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1 Introduction and summary 

The framework of Libor interest rate modeling, initially developed by [T5] , [7] , 
and [13] almost two decades ago, is still considered to be the universal tool for 
evaluation of structured interest rate products. One of the main reasons for this 
is the great flexibility of the Libor framework: It allows to include many sources 
of randomness of different type, such as Brownian motions, Levy processes, or 
even more general semimartingales (see e.g. |15j). Subsequently, these random 
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sources may be connected with different types of volatility structures, such as 
stochastic volatility, local volatility, or deterministic volatilities. In spite of this 
flexibility, the design of a Libor model that can be calibrated in a feasible way 
to a (in some sense) complete set of liquid market quotes (e.g. caps and swap- 
tions for different strikes and different maturities), remains a delicate problem 
however. In its early version, the Libor model was usually driven by a set of 
Brownian motions and equipped with some deterministic volatility structure. 
These Libor models, termed market models, where quite popular because they 
allow for analytic cap(let) pricing and (approximate) analytic swaption pricing 
via Black 76 type formulas. However, a main drawback of these Libor market 
models is that they cannot match implied volatility "smile/skew" behavior ob- 
served in the cap and swap markets. Moreover, these smile/skew effects became 
ever more pronounced over the years. 

For incorporating smile/skew behavior into the Libor model several propos- 
als have been made, for example, the Constant Elasticity of Variance (CEV) 
based extension of the Libor market model by [T], and the displaced diffusion 
Libor market model by [XS] . The implied volatility patterns produced by these 
two approaches have the problem that they are of monotonic nature, so only 
positive or negative skew effects can be imaged. Brigo and Mercurio propose 
in [B] a local volatility model consistent with a mixture of lognormal transition 
densities and some variations on this. One of the problems in this approach is 
the rather complicated volatility structure necessary for Monte Carlo simulation 
of the model in some fixed (e.g. terminal) measure, and the limited flexibility for 
matching too pronounced smile/skew market data. One further line of research 
on smile/skew explaining Libor models concentrates on Libor models driven by 
compound Poisson processes or even infinite activity Levy processes [S]. 
Particularly, in [1] a specifically structured jump driven Libor model is devel- 
oped that allows for feasible sequential calibration to cap volatility-strike data 
for a whole system of maturities. Generally speaking, however, Monte Carlo 
simulation of jump driven Libor models is rather troublesome and expensive 
due to an unavoidable complicated drift term. Recently, in |21j an improve- 
ment is established in this respect, by constructing Levy approximations to this 
Libor drift. In the work of [23J a Heston version of the Libor market model 
is proposed. In the dynamics of this model, which is related to the models in 
[20] and [2] , the volatility of each forward Libor Li (spanning over the interval 
[Tj,T i+1 ]) contains a common stochastic volatility factor y/v where v is a Cox- 
Ingersoll-Ross type square-root process, correlated with Libor driving Brownian 
motions. Moreover, [23. shows that their model has strong potential to produce 
smiles and skews (in particular due to the correlation of v), and they present 
Fourier based quasi analytic approximation methods for the pricing of caps and 
swaptions. Therefore, in some sense the paper of [23] may be considered as a 
first important step towards stochastic volatility Libor modeling. Nonetheless, 
what is missing in this article and in most of the works mentioned above is 
the assessment of the capability of the respective models to be calibrated to 
a larger system of market quotes, including the cap(let) volatility-strike (short 
capvola-strike) panels for a whole system of maturities. In particular, it turned 
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out that only one common volatility factor as in the model of [23) may not 
be sufficient for matching a larger set of cap volatility-strike panels that vary 
significantly over different maturities. The reason is clear: A single stochastic 
volatility factor determines a specific volatility-strike profile that may be consis- 
tent with the market profile over one or some more maturities, but may not over 
a complete tenor structure spanning twenty years for example. As a way out, [2| 
designed a multi-factor stochastic volatility model involving a Brownian motion 
W := (Wk)i<k<n, where the dimension of W is equal to the number of Libors, 
and each component is weighted with a (generally different) square-root type 
stochastic factor v^, and deterministic loading factor (3^, leading to a stochastic 
structure 

dL n_1 

— - = ...dt + J2&kV^dW k , l<i<n, (f) 

h— i 

for the forward Libor Li (under some particular measure). The technical advan- 
tage of this approach is that, after standard freezing of the respective Libors 
Lj to Lj(0) in the drift of the dynamics of L{, a pure affine Libor dynamics 
is produced, and as a consequence, caps and swaptions can be priced quasi- 
analytically by a straightforward extension of the pricing methods in [23] . On 
the other hand, the model of [3] allows for much greater flexibility with regard 
to calibration to a full system of capvola-strike-maturity data. Of course the 
latter doesn't come as a surprise since |lj is in fact a generalization of the 
model in [23] (that is retrieved by taking v k = v). Essentially, in [3] the volatil- 
ity processes in ([!]) are calibrated sequentially to the capvola-strike data in the 
following way. One calibrates the process to the (last) vola-strike panel 

due to L n -\. Next one calibrates v n -2 to the vola-strike data involving L„_2 
with already being identified, and one so works all the way back. After 

carrying out many calibration tests with the model in [3] it turned out that 
the calibration works well as long as there are no big structural movements in 
the capvola-strike patterns when going one step down from 2} to Indeed, 
in the particular case when a larger number of volatility processes are already 
identified, say v^, ...,W4o with n = 41 for an instance, then a single additional 
volatility process V4 may not be able to match a panel at T4 with a sudden 
strongly deviating vola-strike profile. In fact, such breaks in the vola-strike 
patterns where quite typical during the crisis. In this paper we present a new 
flexible multi-factor stochastic volatility Libor model that resolves this problem 
and remains robust even in more critical financial times. 

The central theme of the present paper is a generalization of the Wu-Zhang 
model in the following direction, i.e. we study processes 

^ = ...dt + ^ViPjdW, l<i<n, (2) 
Li 

where by taking Vi = v the Wu-Zhang model is retrieved again. In contrast 
to the structure ([!]), the danger of cumulative cementation of the model in a 
backward recursive calibration is abandoned. Moreover, the dimension of W 
is not strongly restricted anymore to the number of Libors, in order to render 
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a recursive calibration as in Q. However, several technical issues have to be 
resolved. As a main point, even after standard Libor freezing in the drift of the 
full stochastic differential equation (SDE) corresponding to Q, we do not have 
an affine Libor model as in [23] and [3J anymore. That is, the Fourier based 
quasi- analytical approximation for caps doesn't carry over directly. The same 
complication shows up when one attempts to derive an approximate affine swap 
market model from ([2| in order to derive quasi-analytical (Fourier based) swap- 
tion approximations. As a solution we will nevertheless construct affine Libor 
approximations to ([2| and affine swap rate approximations connected with ([2]), 
that allow for quasi-analytical cap and swaption pricing again. But, the price 
we have to pay is that these approximations are typically (a bit) less accurate 
than the ones in the setting of [33J and [3J . Careful tests reveal that the approxi- 
mation procedures developed in this paper are accurate enough for our purposes 
however. The bottom line and justification of our new approach is the following 
"philiosophical" point of view. 

A modeling package that contains only moderately accurate procedures for cal- 
ibrating to liquid market quotes (e.g. accuracy ~ 1%J. but, which is able to 
achieve an adequate fitting error (e.g. ~ 3% due to the 1% off pricing methods) 
in an efficient way, is highly preferable in comparison to a modeling package that 
contains very accurate pricing procedures for calibration (e.g. < 0.2% accurate), 
but, which is unable to achieve an adequate fitting error (e.g. ~ 10%,) despite 
of the accurate pricing formulas. 

Indeed, the former package achieves implicitly a fitting quality with respect 
to the "true model" of about 4%, while the latter package remains left at an 
unsatisfactory fit of ~ 10.2%. Further, for completeness, we extend the structure 
([2| with a standard Gaussian part and with displacement factors like in |16j . 
and consider the structure 



where now W and W are independent standard Brownian motions, ji are deter- 
ministic factor loadings and <Xi are displacement constants for 1 < i < n. From 
a technical point of view this extension goes through without any difficulties, 
neither with regard to the approximate pricing formulas, nor with regard to 
the new calibration procedure. From a practical point of view it enlarges the 
flexibility of the model, but in any particular case, the user can follow her taste 
and may set 7, = 0, or cti = 0, or both. 

As a final introductory note we underline that the cap and swaption ap- 
proximation procedures proposed in this paper can be performed by (inverse) 
Fast Fourier Transformation (FFT), and are thus rather fast. However, as an 
alternative, the recently developed closed form approximation for put/call op- 
tions in a Hcston model from [5], may be straightforwardly adapted to closed 
form cap(let) and swaption pricing formulas in the context of the (approximate) 
affine stochastic volatility Libor and swap rate model here presented. Although 




1 < i < n, 



(3) 
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we consider a detailed treatment here beyond scope, we anticipate that the 
present stochastic volatility Libor model equipped with these formulas might 
be considered an alternative to so called S ABR Libor models (cf. [19] , [12] and 
the references therein). While SABR based models gain popularity because of 
their closed form approximations for vanilla options based on (small time) heat 
kernel expansions, they are also criticized somehow, for instance, because of 
their typically non mean reverting stochastic volatilities. 

2 Recap of Wiener driven Libor modeling 

Let us fix a sequence of tenor dates =: T < T\ < ... < T n , called a tenor 
structure. For each tenor date we consider a zero bond processes B i: i = 
1, . . . , n, where each Bi lives on the interval [0, Ti] and ends with its face value 
Bi(Ti) = 1. A system of forward Libors on the given tenor structure is now 
defined by 

^"-sGsSr 1 )' »£'£T„l & i<n, (4) 

where the periods 5i := Ti + ± — Ti, i = 1, . . . , n— 1, between two consecutive tenor 
dates are termed day-count fractions. In fact, Li may be seen as the annualized 
effective rate due to a forward rate agreement for the period [Ti, Ti + {\ contracted 
at time t. According to this agreement the interest 6iLi(Ti) on the notional 1 is 
to be settled or payed at T i+1 . 

In this article we consider a framework where the Libor defining zero-bonds 
(-^i)i=i n are adapted processes that live on a filtered probability space (f2, 
(•?t)o<t<T 00 , P), where > T n is some finite time horizon and the filtration 
(J-j) is generated by some d-dimensional standard Brownian motion W. Under 
some further mild technical conditions (see [T3] and [TS] for details) there now 
exists for each i, < i < n, an K^-valued predictable volatility process F^ such 
that the Libor dynamics are given by 

IT 11 — 1 r j 

where W*-™- 1 is an equivalent standard Brownian motion under the terminal 
numeraire measure P n induced by the terminal zero coupon bond B n . That is, 
for all j, Bj/B n are P„-martingales. (In this paper we do not dwell on issues 
concerning local versus true martingales.) For some general fixed i, 1 < i < n 
we may consider instead the numeraire measure Pj+i induced by the bond -Bj+i, 
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and then for 1 < j ' < i we obtain from |5]) the dynamics 



dL i = rJ I y 4 f fc r k dt + dww 



= : - E T^VrJr^ + rJdw^ 1 ), i<j<*. (6) 

Since due to Q Lj is a martingale under Pj+i, it automatically follows that 
is a standard Brownian motion under the equivalent measure 
Pi+i. Finally we note that in the case where the Tj are deterministic we have 
the well documented Libor Market Model (LMM) (see for example [5] and [35] 
and the references therein). 



3 A new expiry-wise stochastic volatility model 
with displacement 

The general representation ([5| for the Libor dynamics will now be structured 
towards a multi-factor stochastic volatility model of type (J3|) . Let us take 













H 





w (n) 



dvj = Kj(9j — Vj)dt + y/Vj (aJdW^ 



where 



rJdW {n) 



/(0) = 



J ■■ 



(7) 



where W^ n \ W^ n \ W "' are mutually independent standard Brownian motions 
with dimensions to, to, and, to, respectively, with m + to + to = d. Further, 
for 1 < j < n, f3j and jj are loading factors (in K m and M m respectively) to be 
specified below, and Vj are square-root volatility processes with parameters Kj 
(mean reversion speed), 9j (mean reversion level), and a and tf are deterministic 
"vol of vol" factor loadings (in M. m and R m respectively), where (for convenience) 



(8) 



We thus get 



dLj ^ 5k Lk 

~L~~~ . 1 + S k L_ 

Vj]3jdW {n) +jJdW {n) , 



(ft] Pk^/VjVk + •yJ'Yk) dt 



(9) 
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together with §f§. We next set 



7, = ^7i, & = ^&, (10) 



for deterministic loading factors j3j and 7j (in K m and E m respectively), and 
displacement constants otj, 1 < j < n, and we obtain from ([9]), 

+ JvjPjdWM + 7 7dt? (n) , (11) 

i.e. the new multi-factor stochastic volatility Libor model with displacement and 
stochastic volatilities driven by Q. By applying Ito's formula to the log-Libors, 
(11) becomes 

1 2 1 2 
din (Lj + ctj) = -- Itj-I dt - - Vj dt 

E5 k (L k +a k ) , T T v 

rrr Ll Ui ^ + P 3 h^k) dt 

k=j+l 

+ ^jl3jdW {n) + jJdW^l (12) 



In Section Q we propose a pragmatic approximation that allows for quasi- 
analytical caplet pricing in the context of to (12 1. 

Instantaneous correlations 

For the mutual instantaneous Libor correlations we have 

dLj dLji -p ~ T ~ 

Cov L . L ■- 



7j T 7i' + y/VjVrPjPi' 



V\v\ 2 + Vi\0i\Whr\ 2 + Vi>\Pi>\ 2 ' 

which yields for 7 = 0, Corz^, = and for /3 = 0, Cor^,^, = j^jp^j 

as usual. For the instantaneous correlations between Libors and the stochastic 
volatilities we have 



" ' ' (13) 
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For 7 = we thus obtain 

Cor L . v = 



For the mutual instantaneous correlations between the stochastic volatilities we 
get 

Cor„ v ■= dVj ' dVj ' = a J Tg ■'^ , + °JlL 

3 '' 3 y/dvj ■ dvj^dVj' ■ dvji £ j £ j' 

3.1 Discussion of the Wu-Zhang model as a special case 

Let us take as a special case 7 = 0, otj = 0, Kj = k, 9j = 9, and for some fixed 
unit vectors e vo i G M m , e vo i £ R m , Oj = epe vo i, , ~&j = e^/l — p 2 e vo i where p is a 
fixed correlation constant, — 1 < p < 1. We now are in the setting of Wu-Zhang 



since all volatility processes coincide, i.e. vy = v, and (13) becomes 

Cor£ i)t y = Gov L] . v = peJe vo i, (14) 



where /3j = \/3j\ej with ej G K m . We note that (14 1 reflects a short coming of 
the Wu-Zhang model. The instantaneous correlations between the Libor Lj and 
the common stochastic volatility factor may not be chosen for each j as freely 
as somehow eqn (2.9) from 23 suggests, and we have | Corj^. ^ | < le^e^o;] in 
particular! From another point of view, for realistic uniform skew behavior one 
needs CorL jiV < for all j, so that ej e vo i has to have at least a fixed sign and 
may not become too small for all j. This in turn implies severe restrictions on 
the mutual Libor correlation structure which is usually taken to be an input. 

As an intermediate extension of the Wu-Zhang model above we may consider 
the case 7 = 0, a.j = 0, and then for some unit vectors e vo i € K m , e vo i £ M. m , 



we take <Xj = Ejpje vo i, , <7j = Sj^jl — p^e vo i where pj are fixed correlation 
constants, — 1 < pj < 1, depending on j, and mean reverting speed and level 
may depend on j also. We then have 

CoT L]iV:j = Pjeje voh 

hence for each particular j, any correlation dominated by |eje„ i| may be at- 
tained. Furthermore, as a main feature of the multi-factor model (|^-([l2"]), we 



may have full flexibility regarding the correlations (14), by the structure given 
in Section l4~3l 



Remark 1 If otj = 0, a Libor market model is retrieved by taking j3j = 0, or 
by taking Vj(0) = 6j = 1, <jj = cfj = 0. A further reason for including the 
LMM term 7^ dW in the Libor noise might be to have some extra freedom for 
calibrating to swaptions due to the fact that caplet prices only depend on \ jj \ . 
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4 Approximate caplet pricing and calibration 



For quasi-analytical caplet pricing we will construct an (approximate) charac- 
teristic function of Lj under Pj+i- Let us write (11 1 as 



k=j+i 



5 k (L k + Ofe) 
1 + S k L k 



+ 7*' 



dW {n) 



n-X 

£ 

k=j+i 



8 k (L k + a k ) 



1 + S k L k 



Ikdt 



Since Lj is a martingale under Pj+x, we necessarily have that dW^ +1 ^ and 
are standard Brownian motions under Pj+x- Since the covariation pro- 
cesses (W {n) ,Bj) = for all j, it follows that dW b+1) = dW (n) for all j (cf. 
[25] and [3]). The dynamics of the stochastic volatility process Vj under Pj+x 
can thus be written as 



dvj — Kj(9j — Vj)dt + ^/vj^j dW 



+V v i a 3 



fe=3'+l 



S k ( L k + a k ) 
+ SkLk 



Kj(9 j -Vj)+ \ 



5 k {L k + a k ) T 

Qj Pk^/vjvk dt 



k=j+x 



+ S k L k 3 



(*) 



Thus, in order to obtain approximate affine dynamics for Vj it is enough to 
approximate (*) with an expression that is affine in Vj. Let us therefore consider 
the pragmatic approximation 



y/VjVk 



v k Evj 



VjEvk 



Ev„ 



E Vj 



(15) 



(note that Ev k = 9 k due to the initial condition i n In the Wu-Zhang 

setting we have Vj = v and thus, strict equality in (15) appears. Combining 
(15) and usual freezing of Libors in (*) then leads to the following approximate 
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volatility dynamics, 
dvj rs KjOjdt + 

With 



n-l 



% + E 

fc=i+i 



5 k (L k + a k ) 



1 + 5 k L k 
v](a]dWU + V +ajdW U+1) ) 



n-l 

E 





S k (L k + a k ) 


M 


1 + 6 k L k 



(16) 



we thus obtain from (12 1 the approximative system 



din (Lj + a,) = -1 | 7j | 2 * - i Wj | 2 di + N /u7/?7<iW^' +1 ) + j] 'dW U+1) , 



d Vj = K <f +1) (e ( f +1) - Vi ) ,// 



v](crjdW U+1) +ajdW 



(17) 
v j (0) = 9 j . 



Now the main point is that, if moreover /3,-, tr,, and cr,- are constant in time 



(piece-wise constant would be enough in fact), (17 1 is an afhne structure that 
allows for Fourier based (approximate) caplet pricing. 

4.1 Caplet pricing via characteristic function 

In general the price of a Tj-caplet with strike K is given by 

Cj(K) = SjBj+^Ej+^LjiTj) - K)+ 

= Bj+^SiEj+xiLjiTi) + a, - (K + aj ))+ 

= B^xmEj+xdLjiO) + a,) e 1 " - (K + «,))+ 



B j+1 (p)5 j E j+1 {L* u *(Q)e 



K 



disp\-\- 



We may thus apply the Carr-Madan Fourier pricing method (outlined in the 
next subsection) for caplets using 

disp ^disp / n \ rs-disp 



where the characteristic function 

<pf+f(z;v) := E j+1 



iz In 



Vj(0) = v 



(18) 
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may be obtained as follows. Let us abbreviate for fixed j, X°- X ' v (t) : = In Lj lsp (t) = 
In (Lj(t) + ctj) with X°' x ' v (0) = In if sp (0) = In (Lj(0) + a 3 ) =: x, and V°' x ' v (t) := 
Vj(t) with V°' x ' v {0) = Vj(0) =: v. Then by ^ (using Js])), the generator of the 
vector process (X, V) is given by 

*:-^-(-iw-i.iw')! + ^frr-.)! 

1 /. |2 ,„ l2 \ <9 2 Xn a 2 1 2 <9 2 
+ 2(17,1 +« 1^1)^ + ^7/3^ + 5^. 

Let p (z, z' ; 4, x, u) satisfy the Cauchy initial value problem 

f t =Ap, p{z,z'-^x,v)=e< zx+z ' v ). (19) 

Then 

p(z,z' ;t,x,v) = Ee i (* x °^ v W +z ' v °'" v W). 

We are only interested in the solution for z' = 0. Let us therefore consider the 
ansatz 

p (z ; t, x, v) = exp (A(z] t) + Bq(z; t)x + B(z; t)v) 

with 

Af>;0)=0, B a (z;0) = iz, B(z;0)=0. (20) 
Substitution in |l9|) yields, 



dA dBo BB \ ( 1 a 1 2 \ 

+ K (, + D (^-,)B + i(| 7j | 2 + ,|/3 J | 2 )Bo 2 
+ vaJl3 j B B+~s 2 j vB 2 , 
and we get the Riccati system 

St 
5B 



. w ~|7i| a Bb + «? +1) ^ +1) 5 + i|7i| a Bg 



* ° 



Taking into account ( 20 ) we get 
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It is well known (see [13]) that this system can be explicitly solved, but de- 
pending on the chosen branch of the complex logarithm one may have different 
representations for its solution. We follow Lord and KahPs representation due 
to the principal branch, see [17j^ and obtain 

a,- + di 1 — e djt 

B (z;t) = 3 I J - -r- t 



and 

A(z;t) 

with 



— — 12 + Z 



ds+ 



S °3 



(cbj — dj) t — 2 In 



9j 



1 ~9j 



_ .0+1) 



izaj (3j 



d 3 = v / a 2 + |/3,| 2 (iz + z2) s 2 
9j 



"... - d j 



Resuming, by taking t = Tj we get for (18), 



= exp(A(z;T j ) + B(z;T j )v y ) exp (iz + z 2 ) jf 



\lj\ 2 ds 



(21) 



with 



B(z;Tj) 



A(z;t) :-- 



ij+dj l-e d i T i 
£j 1 — gje d i T i 



and 



(dj -dj)Tj -2 In 



— r\ T 

e 33 — gj 



Carr Sc Madan inversion formula 

Following Carr and Madan [5], the Tj-caplet price is now obtained by the in- 
version formula, 



C 3 (K) = 5 J -B i+1 (0)(if p (0) - iff p )++ 



2tt 



2 



z(z-i) 



dz, (22) 



In a personal communication, Roger Lord confirmed a typo in the published version and 
so referred to the preprint version. 
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where is given by (21 1 and we recall that Vj(0) = Oj. The integrand in 

(22) decays with order z~ z if \z\ —> oo, which is rather slow from a numerical 
point of view. It is therefore advantageous to modify the inver sion formula in 
the following way. Let (fAl sp be the characteristic function (18) due to the 
Black model, 

hf sp {T 3 ) = L<f s P(0)e-H- B ) 2 ^+- B V^, ? e N(0, 1) 

in the measure Pj+i, with a certain suitably chosen volatility . We then have 
(cf. Black's 76 formula) 



E 



(if p (T,) 



K 



disp 



B(Lf sp {0),T J ,a B ,K dlsp ) 7 



where 



B(L, T, a, K) := LM (d + ) — KM (gL) , with 

In4± \a 2 T 

d± := — - — t4 , and 

ay/T 

= e -K- B ) 2 ^(- 2 +-). 

Now applying Carr and Madan's formula to the Black model yields 
Cf{K) := 5 3 B 3+1 (0)B(Lf sp (0),T 3 ,a B ,Kf sp ) = 6 3 B 3+1 (0)(Lf sp (0) - Kf sp )+ 

(23) 



2tt 



z(z — i) 



— iz In 



] dz, 



and by subtracting ( 23 ) from ( 22 ) we get 



C j (K) = Cf(K)+ 



(24) 



2tt 



oq B,disp 



disp i 



z(z — i) 



The latter inversion formula is usually much more efficient since typically the 
integrand decays much faster than in ( 22 ) . 



4.2 Putting the caplet approximation to the test 



We now test the accuracy of the Fourier based caplet pricing method ( 24 ) via the 
approximative characteristic function (21). In this respect we compare, for each 
particular j, the simulation price of the "true" model (11) with the simulation 
price due to the model obtained by replacing each volatility dynamics Vk, k =/= j, 
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with the process Vj, yielding a Wu-Zhang type approximation depending on j 
in fact. In turn, the Fourier based Tj-caplet price approximation is known to 
be a very accurate approximation to the j-linked Wu-Zhang model, as already 
documented in |23j . 

The initial Libor rates are stripped from a given spot interest rate curve (see 
Table [I]). In the test model we drop the Gaussian part, i.e. jj = 0, and also 
assume that no displacement is in force, i.e. a = 0. We choose Sj = T}+i — Tj 



1 .() and we put Q and (11 1 according to Section 4.3 where 
f3j = 0.15ej, such that 



6„- 6-j 



_ -0.073|T 4 - 



-Tj\ 



(25) 



and the other parameters are given in Table [Tj The orthonormal vectors ej 
are obtained by a Cholesky decomposition of the correlation matrix (r^). The 
parameters for the stochastic volatility processes are taken to be representative 
for a typical calibration. In particular they are chosen in such a way that the 
Feller condition 2nd > a 2 is violated. The mean reversion levels are uniformly 
set to 6j = 1. We compare caplet prices due to the "true" model and the 
approximative one, by Monte Carlo simulation based on 30, 000 simulated paths 
(Table [2). 



j 


Pi 


Kj 


£ 3 


Bj(0) 


Lj(0) 


1 


-0 


70 


4 


00000000 


3.00000000 


0.971717 


0.0332468 


2 


-0 


70 


3 


95918367 


2.97959184 


0.94045 


0.0257067 


3 


-0 


70 


3 


91836735 


2.95918367 


0.91688 


0.0195338 


4 


-0 


70 


3 


87755102 


2.93877551 


0.899313 


0.0235296 


5 


-0 


70 


3 


83673469 


2.91836735 


0.878639 


0.0278511 


6 


-0 


70 


3 


79591837 


2.89795918 


0.854831 


0.0258653 


7 


-0 


70 


3 


75510204 


2.87755102 


0.833278 


0.02359 


8 


-0 


70 


3 


71428571 


2.85714286 


0.814074 


0.0237439 


9 


-0 


70 


3 


67346939 


2.83673469 


0.795193 


0.0240497 


10 


-0 


70 


3 


63265306 


2.81632653 


0.776518 


0.023694 


11 


-0 


70 


3 


59183673 


2.79591837 


0.758545 


0.0234799 


12 


-0 


70 


3 


55102041 


2.77551020 


0.741143 


0.0236513 


13 


-0 


70 


3 


51020408 


2.75510204 


0.724019 


0.0238636 


14 


-0 


70 


3 


46938776 


2.73469388 


0.707144 


0.0240064 


15 


-0 


70 


3 


42857143 


2.71428571 


0.690566 


0.0241881 


16 


-0 


70 


3 


38775510 


2.69387755 
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0.0244311 


17 


-0 


70 


3 


34693878 


2.67346939 


0.658177 


0.0246647 


18 


-0 


70 


3 


30612245 


2.65306122 


0.642334 


0.024855 


19 


-0 


70 


3 


26530612 


2.63265306 


0.626756 


0.0249485 



Table 1: Parameters of the Libor model, present values and initial Libor rates, 
terminal bond i?2o(0) = 0.6115. 



The numerical results show that (17 1 approximates very accurately the true 
model dynamics Q and (11). Indeed, the absolute price deviations are of 
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magnitudes within basis points, with a well behaved relative error for ITM 
(in-the-money) and ATM (at-the-money) contracts. The relative errors become 
somewhat larger for OTM contracts, but OTM (out-of-the-money) caplet prices 
are typically very low (close to worthlessness) so that relative errors stemming 



from approximation (151, ( |17[ ) are intrinsically unstable (for any "good" approx- 
imation in fact). 
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(9.43e-05) 
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0.0430 




0.020 


0.0055 


(8.62e-05) 


0.0051 


(8.08e-05) 


3.25e-04 


0.0509 




0.025 


0.0040 


(7.72e-05) 


0.0037 


(7.17e-05) 


3.12e-04 


0.0773 




0.030 


0.0029 


(6.81e-05) 


0.0026 


(6.26e-05) 


2.84e-04 


0.0963 



Table 2: Simulation results for caplets. 



4.3 Further structuring and calibration 

As part of the model, we choose a fixed LMM part jj of the Libor structure. 
This part may be obtained from an LMM calibration, eventually weighted with 
some factor for instance or, if enough flexibility is left for our purposes, we may 
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set 7j = 0. The loadings f3j are also assumed to be chosen in advance. We 
further take m = 1 in Q, and for pj, —1 < pj < 1, we take o~j —: EjPjej, 

where j3j =: |/3j|ej, and so Oj =: «/l — p 2 £j. Note that in principle we have no 



restrictions on pj conferred to the Wu-Zhang case (see Section 3.1). Then (13) 
becomes 

( '"'•/. ••• • = P3' e l e f = Pf r n' 
with Tjji := sjej', and in particular we have Cor^. ^ = pj. For the mutual 
correlations between the volatility processes we so have 



Cor„ . >ty = Pjp r r, jr + yjl-p* \J l ~ p} ■ (26) 

In any case the scalars Kj , 6j , pj , Sj , and the loadings have to be time indepen- 
dent, in order to invoke standard square-root volatility processes. In principle 
piecewise constant t j3j (t) will allow for Fourier based caplet pricing later on, 
but for simplicity we assume henceforth that the /3j are also time independent. 

Remark 2 In practice it turns out that the pj are negative overall in order to 
produce a skew. Let us assume for simplicity that we could fit the data with a 



uniform (negative or positive) p. Then (26) implies Cor Vj _ Vjl = 1 — p (1 — rjj 



> 1 — p 2 , assuming that mutual Libor correlations rjj> are nonnegative. This 
means that mutual correlations between volatility processes are typically high 
(> 0.5 for p — 0.7), and even close to 1 when j' is close to j. 

4.4 Calibration to caplet volatility-strike-maturity 

We will now illustrate a typical calibration test of the stochastic volatility Libor 
model in its terminal measure to market cap-strike data. The test is carried 
out for EurlBOR market data from September 20, 2010, based on a twenty year 
semi-annual tenor structure. For simplicity, the displacements and the Gaussian 
part where taken to be zero, i.e. ay = 0, 7,- = 0, and as further input parameters 
we took 6i = 1, and from a Cholesky decomposition according to ej ej = r.^ 
_ g-o.ii8lTi-T.fl_ p or gach maturity Tj, the parameters 

\f3j\ , Kj, Ej, Pj, 

where next calibrated to the caplet price-strike panel corresponding to Tj, ob- 
tained from the market data. This calibration involves a minimum search of a 
standard averaged relative error functional based on the FFT pricing formula 
( 22 ) due to the characteristic function ( 21 ). Each trial Kj (which is restricted to 



with pj , are subsequently plugged into ( 2T ) 



> 0) induces a Kj 1 '* and 6j via (16) (recall that Q{ = 1) which, together 

The implied volatility patterns due 



to the calibration as well as the calibrated parameters are depicted in Figure [TJ 
Concluding we may say that we obtained a satisfactory model fit with robustly 
behaving parameters when moving from one maturity to the other. Optically 
the fits for small strikes, hence deep ITM caplets may look a little bit off over- 
all. However, this is only appearance because our algorithm calibrates to caplet 
prices, while implied volatilities are badly conditioned for deep ITM strikes. 
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Figure 1: Implied caplet volatilities due to market data vs. calibrated model 

17 



5 Swap rate dynamics and approximate swap- 
tion pricing 

5.1 Swap contracts and dynamics under swap measures 

An interest rate swap is a contract to exchange a series of floating interest 
payments in return for a series of fixed rate payments. Consider a series of 
payment dates between T p+ i and T q , q > p. At each time j = p, . . ■ , q — 1, 

the fixed leg of a (standard) swap pays 5jK, whereas in return the floating leg 
pays SjLj(Tj) with Lj(Tj) being the spot Libor rate. Consequently, the time 
t— value of the interest rate swap (with t < T p ) is 

9-1 

J2^B j+1 (t)(L^(t)-K). 
j=p 

The swap rate S Ptq (t) is defined to be the value of K for which the present value 
of the contract is zero. We thus have 



E?4 8 3 B 3+l (t)L 3 {t) _ Bp (t) - B q {t) 



9-1 



(27) 



So S Pj q is a martingale under the probability measure P p ,q, induced by the 
annuity numeraire 

9-1 

BpAt) :=X>£j+i(*). 

3=P 



From (11) it follows that 

dS p>q (t) = S p . q (t)Al q (t)dW^(t), 



(28) 



where :— {W p,q , W p,q ) is standard Brownian motion under P p , q , and 

where 



j=p J \ i=j 



B„ 



Bp ~ Bg 



lj 



SjBi+i 



(29) 

The derivation hereof is given in Appendix[6j We further have (see Appendix[6]), 

(30) 



9-1 n-1 x IT \ \ 

dW™ = dWW - dt £ < 9 £ ^ fe ~ r ' ! ' 1 



7fc 



By ( 28 ) we thus get 



d In 5, 



P,9 



P,9 

1 



5, 



(31) 



/J, 9 



|A p , g | 2 ^ + A^ 9 dW^, 
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where by ( 29 ) we may write 



with 



9-1 
i=p 



1 + SjLj 



7j 



'9-1 



S,, 



Bn — B„ 



(32) 



J=3 



B„ 



3i 

B n 



(Cf. [22], (1.35), and (1.38) so we have that £?' 9 (0) wf' 9 (0) with equality 
when the yield curve is flat; hence the are approximate weights also.) 



5.2 Approximate affine swap rate dynamics 

In order to approximate the swap rate process with a pure square-root volatility 
process we introduce the process 

dv v,q = K P,Q(0P,Q _ v p,i) dt + ^f-p^ ( a J q dW {n) + aJ tq dW (n) ) , v p - q {0) = 6> p <« 

(33) 

with 

9-1 

.= ^^.9(0)5, 
Z=p 
9-1 

:=^ w f^(0) Ki 

9-1 

a M :=£< 9 (0)^. (34) 
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By replacing in (32) all volatility processes Vj with the, in a sense, averaged 
process v p,q , and freezing Libors we arrive at the approximation 



A 



p,i 



q-l 9-1 

EE 

j= p j= p L 7j 



(0) 



where 



Pp,<3 



7p,g 



9-1 

Eft 

3=P 

E-* 

j=p 



(0) and 
(0) 



(note that X)j=p 9 (-^3 /^p-? ~ 1)> hence yielding affine approximative 

swap rate dynamics 

d\nS p , q = ~v pq \P m \ 2 dt~^ \ lp , q \ 2 dt + Vv^Pj tq dW p < q + ll q dW™. (35) 

For the (approximate) dynamics of v p ' q under the annuity Brownian motions 
we replace in (30 1 the processes Vj by their average v p,q , and freeze the Libors 
as usual. From ( 33 ) we then obtain (as in Section |4j it follows again that 

w in) =w p - q ), 



v p - q )dt + Vv p ' q a p \ q dW 

q—1 n—1 r 



+ Vv p ' q <7, 



p.q 



dW p ' q 



/ ^ rf *E E 

;=p k=i+i 



M 

1 1 + 6 k L k 



(0)/3 fc 



By setting 



q-l 

K p ' q := n p ' q - 

i= P 

K p - q Q p - q 



n-1 



p,q fc I 

,7 f e -\ 



k=l+l 



4 ( £fc + a k ) 
+ SkL k 



(0)viqPk 



•,q _ 



we thus have (in approximation) 



dv p ' q = K 



p,qfap,q _ ,,p,q 



- v p < 9 )cft + Vv p - q a p \ q dW p - q + Vv p > q a ptq dW* 



(36) 



5.3 Fourier based swaption pricing 

A (payer) swaption over the period [T p ,T q ] is the option to enter at T p into a 
swap over the period [T p , T q ] with strike K. It follows straightforwardly that the 
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value at time t = is given by 

Swpn P , g ( K ) = B Ptq (0)E p , q \(S M (T p ) - K) + 



(37) 



Thus, after determining the characteristic function for In [S Piq (T p ) / S p . q (0)} we 
may price the option by the Carr-Madan Fourier inversion method, just like we 
did for caplets in Section |4~T) Recalling the analysis from Section |4~1~1 it follows 
immediately that this characteristic function is given by 



<P P ,q(z;v) '■= E p-q 



v p , q (0) = v 



(38) 



exp (A p<q (z;T p ) + B p ^ q (z;T p )v) exp I -- (iz + z 2 ) J \"/ p<q \ ds 



where 



a 4- d 1 — pd p qTp 
P -q\ i p) 2 i d T 

1 — g p _ q e a p-i 1 p 



and 



p,q 

-^p,q0P,q ( 
T2 7^1 T2 1 \®PiQ ^P,q) ^P ^ ^ 



'p,q\ 



'p,q\ 



e «p.^p -g pq 



1 - 9p,< 



with 



^p.q K^'^ \Z(T pq /3 p . q 



dp,q = \j a l,q + l/3p,?| 2 + Z 2 ) (Wp,q\ 2 + Wp,q\ 2 ) 



9p,q 



^p,q ~t~ ^p,q 
^p^q dp 7 q 



Based on (38 1 the (approximate) price of a swaption with maturity T p and 
swaption leg [T p ,T q ] is given by 



Swpn p . q (K) = B Piq {0)E p , q [(S p , q (T p ) - K) + 
« S«;pnf i ,(A-)+ 

^,9(0)6^,5(0) f°° <p® tq (z-i;T p ,6 Piq ) - (p Ptq (z-i;T p ,6 p , q ) _ izln 



(39) 



2- 



z(z — i) 



s p,«j<°) (iz 



In ( 39 ) , is the characteristic function of a corresponding Black model, 

S p , 9 (T p ) = S p . q {Q)e-^<^ T - + <^\ ? G tf(0, 1), 
where up 3 g is a suitably chosen volatility, and 
S^n^if) = B p , q (fl)E M (S p , q (T p ) - K)+ = B M (0)B(S p , q (0),T p , a? K), 



is given by Black's formula (cf. (24)) 
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5.4 Putting the swaption approximation to the test 



In the same spirit as we have tested the caplet price approximation in Section |4~2"1 
we now test the above Fourier based swaption pricing method. For each pair 
(p, q), 1 < p < q < n {q ^ p+1), we replace all volatility processes Vj, p < j < q, 



with Vp^ q given by (33), (34) to obtain in fact a Wu-Zhang related swaption 



approximation model linked to this pair (p, q) . We then compare the simulated 
(p, q , )-swaption price due to the "true" model ( 11 1 and the model with common 
stochastic volatility process (33). In turn, the latter price can be accurately 
approximated by ( 39 ) as shown in [24] . We base the numerical experiments on 



the same data set as in Section 14.21 

In detail, this means that for putting up the "true" and the approximate 
Libor model, the initial Libors are stripped from a given spot rate curve and 
their values are given in Table [l] the Gaussian 7-part is deactivated by putting 
7j = and no displacement is in force by choosing ay = 0. Moreover, the 



parametrization of the correlation structure from Section |4.3| is given by 

r i:j = exp ( - 0.0553IT, - T 3 \) = ejej, (3j = 0.15e,-, 

with the orthonormal vectors resulting from a Cholesky decomposition of 
(rij) and 5j — Tj+± — Tj = 1.0 and 8j = 1 remain valid. All other simulation 
parameters, in particular the Pj's, Kj's and e^'s can be found in Tableland we 
retain the diffusion coefficients 



Pj £ 3 



p]e 3 . 



To gear towards the approximate Libor model, we perform the calculation 
of the weighted volatility parameters n p,q , 9 p ' q , a p,q and a p,q according to (34), 
where the frozen weights wf' q (0) are given in (29), so that the averaged ap- 
proximate volatility process v p,q from ( 33 ) can be simulated. This averaged 



i.e. v 



stochastic volatility is then reinserted into the Libor dynamics ^ , 
tually replaces each expiry- wise volatility Vj,j=p,...,q—1. The simulations 
are carried out using 30, 000 Monte Carlo paths. 

We calculate "true" and approximate swaption prices for the payer swaption 
depicted in (37) for various strike levels and swap legs [T p ,T q ]. The results of 
our numerical experiments are depicted in Table [3] 

The simulation results show that for swaption pricing, the approximate Li- 
bor model under one weighted stochastic volatility v p,q gives a surprisingly good 
fit to the true model dynamics ([9|, Q. Depending on the swap legs, absolute 
price deviations are in the range of basis points (for swaption maturing in two 
and four years) and in the range of ten basis points (for maturity ten years). Re- 
calling that the approximation is somewhat strong as each expiry-wise volatility 
process Vj,j = p, .. . ,q — 1 is replaced by one weighted volatility process v p,q , 
the numerical results reveal however that we get reasonably well behaved ap- 
proximations to the "true" model. 
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(7.1e-05) 


n nn/i n 
U.UU4U 


(o.ye-uo ) 


U.UUUzD 


n n^n 
U.UOU 




0.000 


0.2877 


(4.8e-04) 


0.2866 


(4.8e-04) 


0.00110 


0.003 




0.005 


0.2288 


(4.6e-04) 


0.2277 


(4.6e-04) 


0.00107 


0.004 




0.010 


0.1699 


(4.5e-04) 


0.1689 


(4.4e-04) 


0.00104 


0.006 


r a onl 
[4, 20J 


n m c 
U.Uio 


0.1122 


(4.2e-04) 


U.lllz 


( A f\A\ 

(4.ZC-U4 J 


n nm no 
U.UUlUz 


n nnn 

u.uuy 




0.020 


0.0609 


(3.5e-04) 


0.0600 


(3.5e-04) 


0.00091 


0.015 




0.025 


0.0246 


(2.4e-04) 


0.0241 


(2.4e-04) 


0.00051 


0.020 




0.030 


0.0068 


(1.2e-04) 


0.0067 


(1.2e-04) 


0.00011 


0.016 




0.000 


0.1653 


(4.5e-04) 


0.1638 


(4.4e-04) 


0.00149 


0.009 




0.005 


0.1311 


(4.4e-04) 


0.1297 


(4.3e-04) 


0.00146 


0.011 




0.010 


0.0976 


(4.2e-04) 


0.0961 


(4.1e-04) 


0.00146 


0.015 


[10, 20] 


0.015 


0.0670 


(3.9e-04) 


0.0655 


(3.8e-04) 


0.00147 


0.021 




0.020 


0.0423 


(3.3e-04) 


0.0410 


(3.3e-04) 


0.00137 


0.032 




0.025 


0.0247 


(2.7e-04) 


0.0236 


(2.6e-04) 


0.00137 


0.045 




0.030 


0.0134 


(2.0e-04) 


0.0126 


(1.9e-04) 


0.00081 


0.060 



Table 3: Simulation results for payer swaptions. 



6 Appendix 



The derivation of the swap rate volatility ( 29 1 is essentially given in [55] . But in 
order to match to the present notation and to make reading more convenient, we 
now give a short recap. Let, exclusively in this section, <jj denote the volatility 
of the bond Bj, let fij be the drift of Bj, and A be the market price of risk 
process with respect to the driving Brownian motion VV =(W, W, W). That is, 
in the objective measure the zero bond dynamics are of the form 

dB 

— - = Hjdt + aJdW with (ij = erjA, 
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and where <7j- & = for m + m<k<m + m + m. Following [531 P-17], we may 
write 

9-1 9-1 

dB p,q = E 5 3 dB j+i = ■ ■ + SjB j+ ioJ +1 dW 

j =p j— p 

9-1 

= ---dt + B p , q wf q aj +1 dW. 
j=p 

We thus have by Ito's formula for p < r < q, 



d(B r /B p ^ q ) 

Br I Bp q 



9-1 



■dt+ \aJ-J2wfaJ +1 ]dW 



9-1 



as B r /Bp tq is a P p 9 -martingale. We thus obtain 



dSp : q — d 



J p,i 



Bp 
~B~ 



p.q 



-7-X>r- 



7+i 



B, 
B 



11,11 



9-1 

^ - E < v 7+i 



Bp - E -r « - -#^E < 9 k t - *7 + o 



B 



P,9 



P,9 



dWP' 1 



S. 



p.q 



B„ 



Bp — B q 



9-1 

E 



(-7 - ct 7+0 



B H 



Bp — Bq 



9-1 

E 



< 9 K T - -7m) 



= 5, 



p.q 



9-1 



E «$* (-7 - *J + i) + iAw E < q (-7 - -7) 



9-1 



=: S^A^dW™. 
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Similar to (1.13) in [S3] we get 



9-1 



B 



J=P 
9-1 



B„ 



9-1 



E^E^ _0 >+l) + B ^ B ^{Or-Or+l) 
j—p r—p P ^ r=p 



9-1 



'9-1 



= ^ (°> - CTr+l) E' 



,P,9 



r—p 



-Bp — 



9-1 

E 

r—p L Tr 



V r P r 



1 + 5 r L r 

Further, by (1.27) from [22], it holds that 



'9-1 

E-r 



B p - B q 



9-1 



dW {n) — dW + (A — <7„)eft, and <M p ' q = Xdt wf' q a l+1 dt + dW. 

i=p 



Therefore, we finally have 



9-1 



dW™ = dW {n) + a n dt - J2 wf' q a l+1 dt 



i=r, 



9-1 



i= P 

9—1 n— 1 

= dw [n) +dt j2 < ? E - ^) 

5fe (L fc + ajt) 



z=p fe=i+i 

9—1 n— 1 



fe=;+i 



Ik 
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